clear all

%load success_femot;
load results_famot_MRSA.txt -ascii

save_array = results_famot_MRSA;

save_array(:,3)= save_array(:,3) ./ save_array(:,4);
save_array(:,4) = 0;

summary = sortrows(save_array, [1 2 3])


% format of summary
%
% temperature: [10 15 20 25 30];
%  P_acc : 0.2 0.35 0.5 0.65 0.8 0.95
%  N_md/N_mc : 2/250000 5/100000 16/31250 50/10000 158/3165 500/1000
%  1000/500


% notice that the first two special coordinates are swapped
% this is the peculiar format of isosurface, slice and so on.


Ntemp = 5;
Naccept = 6;
Nratio = 7;

V = zeros(Naccept+1,Ntemp+1,Nratio+1);


% each first index value of TempSlices contain a acceptRatio
% versus NmdNmc image (notice the two last indices have been
% swapped because this is have matlab work with images)

WhichToPlot = 5;

index = 1;
for i = 1 : Ntemp
    temperature(i) = summary(index,1);
      
    for j = 1 : Naccept
        if i == 1
            acceptRatio(j) = summary(index,2);
        end
        
        for k = 1 : Nratio
            if i == 1 & j == 1
                NmdNmc(k) = summary(index,3);
            end
            
            V(j,i,k) = (summary(index,WhichToPlot));
            
            %V(j,i,k) = -log10(summary(index,2))*summary(index,1);
            
            % cutoff to better view step-length
            
%             if V(j,i,k) > 0.05
%                 V(j,i,k) = 0.05;
%             end
            
            % add extra line to circumvent the problem that slice
            % do not slow last values in array
    
            if i == Ntemp & j ~= Naccept & k ~= Nratio
                V(j,i+1,k) = V(j,i,k);
            end
            if i == Ntemp & j == Naccept & k ~= Nratio
                V(j+1,i+1,k) = V(j,i,k);
                V(j+1,i,k) = V(j,i,k);
                V(j,i+1,k) = V(j,i,k);
            end
            if i == Ntemp & j == Naccept & k == Nratio
                V(j+1,i+1,k+1) = V(j,i,k);
                V(j+1,i+1,k) = V(j,i,k);
                V(j+1,i,k+1) = V(j,i,k);
                V(j,i+1,k+1) = V(j,i,k);
                V(j+1,i,k) = V(j,i,k);
                V(j,i+1,k) = V(j,i,k);
                V(j,i,k+1) = V(j,i,k);
                
                temperature(i+1) = temperature(i) + (temperature(i)-temperature(i-1));
                acceptRatio(j+1) = acceptRatio(j) + (acceptRatio(j)-acceptRatio(j-1));
                NmdNmc = log10(NmdNmc);
                NmdNmc(k+1) = NmdNmc(k) + (NmdNmc(k)-NmdNmc(k-1));
            end
            
            if i == Ntemp & j ~= Naccept & k == Nratio
                V(j,i+1,k+1) = V(j,i,k);
                V(j,i,k+1) = V(j,i,k);
                V(j,i+1,k) = V(j,i,k);
            end
            
            if i ~= Ntemp & j == Naccept & k == Nratio
                V(j+1,i,k+1) = V(j,i,k);
                V(j+1,i,k) = V(j,i,k);
                V(j,i,k+1) = V(j,i,k);
            end
            if i ~= Ntemp & j ~= Naccept & k == Nratio
                V(j,i,k+1) = V(j,i,k);
            end
            if i ~= Ntemp & j == Naccept & k ~= Nratio
                V(j+1,i,k) = V(j,i,k);
            end
            
            
            
            % add to TempSlices
            %TempSlices(i,k,j) = summary(index,4);
            
            index = index + 1;
        end
    end
end
    
Ekin_factor = 13/2;

% subplot(1,2,1)
[x,y,z] = meshgrid(Ekin_factor*temperature, acceptRatio, NmdNmc);
H=slice(x,y,z, V, Ekin_factor*[10 15 20 25 30], max(acceptRatio)+1000.05, NmdNmc(1)-300.05)
for i = 1 : length(H)
  set(H(i), 'edgecolor', [0.2 0.2 0.2]);  
end
%slice(x,y,z, V, 500, [0.65 0.8 0.95], NmdNmc(1)-500)
%axis([7 17 0.15 1 -6 2.5])
%slice(x,y,z, V, 11, [0.65 0.8 0.95], NmdNmc(1)-500)
colorbar
%scatter3(temperature, acceptRatio, log10(NmdNmc),V)
%scatter3(x(:),y(:),z(:), V(:))
xlabel('\langle E_{kin} \rangle')
ylabel('p_{a}')
zlabel('log_{10}( N_{md} / N_{mc} )')
axis([Ekin_factor*8 Ekin_factor*32 0.2 1.2 -6 2])
colormap(gray)
hold off
